In [1]:
def matrix_is_square(M):
if len(M) == len(M[0]):
return True
else:
return False
class MatrixIsNotSquareException(Exception):
pass
class MatrixIsSingularException(Exception):
pass
def matrix_id(d):
'''
This function takes an integer d and
returns the d*d identity matrix
matrix_id(3)
[[1, 0, 0], [0, 1, 0], [0, 0, 1]]
'''
return [i*[0]+[1]+[0]*(d-i-1)
for i in range(d)]
def matrix_row_swap(M, upper_row, lower_row):
'''
>>> matrix_row_swap([[1,2],[3,4]],0,1)
[[3, 4], [1, 2]]
'''
M[upper_row], M[lower_row] = M[lower_row], M[upper_row]
return M
def matrix_row_multiplication(M, row_number, scalar):
'''
matrix_row_multiplication([[1,2],[3,4]],0,3)
[[3, 6], [3, 4]]
'''
M[row_number] = [scalar * e for e in M[row_number]]
return M
def matrix_row_addition(M, changed_row, mulitplied_row, scalar):
'''
matrix_row_addition([[1,2],[3,4]], 1, 0, 2)
[[7, 10], [3, 4]]
'''
M[changed_row] = [sum(x) for x in zip(M[changed_row],
[scalar * e for e in M[mulitplied_row]])]
def matrix_inverse(A):
'''
This function calculates the inverse B of a given matirx A
such that AB = BA = I, where I is the identiy matirx of
same dimension as A or B.
>>> matrix_inverse([[2,1],[0,2]])
[[1.0, -0.5], [0.0, 0.5]]
>>> matrix_inverse([[1, 1, 0], \
[0, 1, 1], \
[2, 0, 2]])
[[0.5, -0.5, 0.25], [0.5, 0.5, -0.25], [-0.5, 0.5, 0.25]]
>>> matrix_inverse([[1, 1, 0, 1], \
[0, 1, 1, 1], \
[1, 0, 1, 1], \
[1, 1, 1, 1]])
[[0.0, -1.0, 0.0, 1.0], [0.0, 0.0, -1.0, 1.0], [-1.0, 0.0, 0.0, 1.0], [1.0, 1.0, 1.0, -2.0]]
'''
if not matrix_is_square(A):
raise MatrixIsNotSquareException
B = matrix_id(len(A))
# Make all elemens under the diagonal 0
# For each column in the matrix
for col in range(len(A)):
# Swap the rows such that elements on the diagonal are not zero.
if A[col][col] == 0:
for row in range(col+1, len(A)):
if A[row][col] != 0:
matrix_row_swap(A, col, row)
matrix_row_swap(B, col, row)
break
else:
raise MatrixIsSingularException
for row in range(col+1, len(A)):
if A[row][col] != 0:
scalar = -(A[row][col]/A[col][col])
matrix_row_addition(A, row, col, scalar)
matrix_row_addition(B, row, col, scalar)
# Make all elemens over the diagonal 0 and the diagonal 1
for col in range(len(A)-1, 0, -1):
if A[col][col] != 1:
scalar = 1/A[col][col]
matrix_row_multiplication(A, col, scalar)
matrix_row_multiplication(B, col, scalar)
for row in range(col-1, -1, -1):
if A[row][col] != 0:
scalar = -(A[row][col]/A[col][col])
matrix_row_addition(A, row, col, scalar)
matrix_row_addition(B, row, col, scalar)
return B
if __name__ == '__main__':
import doctest
doctest.testmod()
In [2]:
def fact2(n):
'''
This function takes a positive integer n
and returns its factorial
>>> fact2(0)
1
>>> fact2(1)
1
>>> fact2(8)
40320
'''
assert type(n) == int, \
'd should be an integer it is an '+type(n)
assert n >= 0 , \
'd should be a positive integer it is '+str(n)
if n == 0:
return 1
result = 1
for x in range(1,n+1):
result *= x
return result
def prod(l):
'''
This function takes a list of integers and
returns the product of all elements
>>> prod([1,2,3])
6
>>> prod([-1,2,-3])
6
>>> prod([3,0,2])
0
>>> prod([-3,5.5,1,7])
-115.5
'''
assert type(l) in [list, range], \
'l should be a list or range of numbers'
if len(l) == 1:
return l[0]
if len(l) == 2:
return l[0]*l[1]
else:
mid = len(l)//2
return prod(l[:mid]) * prod(l[mid:])
def fact3(n):
'''
This function takes a positive integer n
and returns its factorial
>>> fact3(0)
1
>>> fact3(1)
1
>>> fact3(8)
40320
'''
assert type(n) == int, \
'd should be an integer it is an '+type(n)
assert n >= 0 , \
'd should be a positive integer it is '+str(n)
if n == 0 or n == 1:
return 1
return prod(range(1,n+1))
if __name__ == '__main__':
import doctest
doctest.testmod()
In [3]:
import time
from functools import wraps
def timeit(func, calls=1):
'''
Decorator to time the execution time of a function
'''
@wraps(func)
def timed(*args):
measurements = []
for call in range(args[-1]):
start = time.perf_counter()
result = func(*args[:-1])
total = time.perf_counter() - start
measurements.append(total)
avg_time = sum(measurements)/calls
return result, total
return timed
In [4]:
import math
timed_fact2 = timeit(fact2)
timed_fact3 = timeit(fact3)
timed_math_fact = timeit(math.factorial)
x = []
time_fact2 = []
time_fact3 = []
time_math_fact = []
for exp in range(1,21):
x.append(2**exp)
time_fact2.append(timed_fact2(2**exp, 10)[1])
time_fact3.append(timed_fact3(2**exp, 10)[1])
time_math_fact.append(timed_math_fact(2**exp, 10)[1])
In [5]:
timed_math_fact(2**20, 10)[1]
Out[5]:
In [7]:
%matplotlib inline
import matplotlib.pyplot as plt
plt.style.use('ggplot')
plt.loglog(x,time_fact2,':', label='fact2', )
plt.loglog(x,time_fact3, '--', label='fact3', )
plt.loglog(x,time_math_fact, '-.', label='math.factorial', )
plt.legend(loc=2)
plt.xlabel('Function input')
plt.ylabel('Execution time [s]')
fig = plt.gcf()
fig.set_size_inches(8,5)
In [14]:
skeleton = '{:^10.5} {:^10.5} {:^10.5}'
print(skeleton.format('frac2', 'frac3', 'mfrac'))
for x in zip(time_fact2, time_fact3, time_math_fact):
print(skeleton.format(*x))